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Abstract 

We determine to order a 3 in the quenched approximation the so-called residual 
mass in the lattice regularisation of the Heavy Quark Effective Theory. We follow 
a gauge-invariant strategy which exploits the fact that this mass term dominates 
the exponential decrease of perturbative Wilson loops at large perimeters. Our 
computational tool is Numerical Stochastic Perturbation Theory. The new coef- 
ficient we compute is crucial to improve the determination of the (MS) mass of 
the 6-quark from lattice simulations of the Heavy Quark Effective Theory. 
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1 Introduction 



Heavy Quark Effective Theory (HQET) is proving to be a particularly useful tool for phe- 
nomenological studies in heavy flavour physics (see [Q for a review and references therein 
to original papers). In this approach physical quantities are expanded as series in inverse 
powers of the heavy quark masses mg >> Aqcd- The Lagrangian density of the effective 
theory may be written as 

£hqet = hD 4 h, (1) 

where D4 is the (gauge) covariant time derivative; h is the bi-spinor describing the only 
degrees of freedom of the heavy quark that are relevant in this approximation. Soon after its 
introduction it was realized that while HQET could appear conceptually simple, there was 
a variety of theoretical subtleties one should face in order to profit from its application. Not 
surprisingly, many of these subtleties have to do with the definition of the heavy quark mass 
rriQ itself. First of all it was realized that the theory must be invariant by a redefinition 
of rriQ which leaves the momentum Pq fixed: 

P» = m Q v» + fc". (2) 

In fact the mass of a quark is not a well defined physical quantity, and the separation of 
Pq into a large and a small (k^ ~ Aqcd) component is to some extent arbitrary. The 
requirement of invariance by redefinition of rriQ leads to the introduction of a "covariant" 
derivative in terms of which the HQET Lagrangian can acquire a mass term for the static 
quark. It was noted however that even so the expansion in 1 / rriQ is not well defined and other 
ambiguities show up. In dimensional regularisation they appear as UV renormalons in the 
matrix elements of HQET and IR renormalons in the coefficient functions that match to QCD 
0. In the lattice regularisation problems show up as non-perturbative power divergences 
in the inverse lattice spacing [|J. The ambiguities in the coefficient functions are expected 
to cancel with those in the HQET matrix elements. In order to perform such a subtraction 
and to obtain a properly defined expansion the expansion parameter itself {ttlq), first of all, 
must be defined in a way free of ambiguities. The pole mass m pole proved not to be a good 
definition [||]. The authors of |5|] introduced a lattice non-perturbative definition for the 
quark mass (the so-called subtracted pole mass) which is free of renormalon ambiguities and 
power divergences 

tjiq = m H — A. (3) 

In the previous formula m# is the mass of an hadron containing the heavy quark, and A is 
the non-perturbative definition of a subtracted binding energy 

A = £ - 8m. (4) 

£ is the "naive" binding energy, computed as the decay constant of the correlation function 
of two axial currents (A4 = h(x)'j4'y 5 q(x), q(x) being a light quark) evaluated with the action 
in Eq. ©: 

Cit) = ^2(A 4 (x,t)A 4 (0,0y) (5) 
~ Z 2 exp(-£t). (6) 
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The other term in Eq. (Q) is a non-perturbative definition of the so-called residual mass 
counterterm. The latter is just the quantity we are concerned with in this paper. It is 
a linearly divergent mass counterterm which quantum corrections generate in the lattice 
regularisation, contrary to dimensional regularisation. In the non perturbative approach we 
have just referred to, it is fixed by the large time behaviour of the quark propagator S%(x, y) 



,j colour indexes): 



5m 



1 i r 1 1 — log 

t^oo a 
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(7) 
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Both mq and A are non perturbative, well defined quantities and may be evaluated by 
numerical simulations on the lattice || [Tj. This is not at all a trivial observation. Both £ 
and 5m are linearly divergent quantities on their own and only the result of the subtraction 
in Eq. (||) is a well defined, divergences free quantity, which also makes well defined Eq. (|3|). 
However, the subtracted pole mass is still a long distance quantity and in order to be able to 
extract physical predictions from lattice HQET it is necessary to match the definition of mn 
to another mass definition such as the mass tTTq of MS. In fact the coefficient functions of the 
1 / mg expansion are typically computed in perturbation theory in dimensional regularisation 
||. The relation between mg and the MS mass of the 6-quark m^(m b ) has been computed 
in H [7] and exploited to compute the mass of the 6-quark. The strategy goes as follows. 
The pole mass (we write m v h ole having in mind the pole mass for the 6-quark) can be related 
to the MS mass in perturbation theory 
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Note that in this relation coefficients D n with n < 2 are known [§, [|. By matching the 
propagator in QCD to its lattice HQET counterpart (we closely follow 0) one can then 
relate the pole mass to the binding energy and the physical mass of an hadron, for example 
the B meson mass (Mb)- This relation is just the perturbative version of Eq. (|3]) 
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0(l/m b ) = M B -£ + 5m + 0(l/m b ) 



(9) 



This time the residual mass term (5m) has naturally emerged from the perturbative compu- 
tation of the heavy quark propagator on the lattice. We slightly change notation to stress 
that 5m is the perturbative residual mass. By substitution one gets 
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(10) 



Note that the last expression asks for both a non-perturbative (£) and a perturbative (5m) 
computation. 5m is still in charge of cancelling the linear divergence of £, even if now the 
cancellation is in perturbation theory, so that it is incomplete. What one also gets from 
Eq. ( |10D is a delicate cancellation of renormalon ambiguities as well: the renormalon in the 
expansion of 5m cancels the one in the perturbative relation between the pole and MS masses 
and this is the key issue. This cancellation asks first of all for the use of the same coupling 
in both perturbative expansions. It also asks for high orders in the perturbative expansions. 
Since, as we have already said, D 2 is known, this means that X 2 has to be computed. In 
the end it turns out that the ignorance of this term is the main source of error in the 
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final determination of m^{rn^). This observation ends this introductory review which sets 
the motivations for our computation: the purpose of this paper is to compute 5m to 
order (X 2 ) in the quenched approximation. Our computational tool is Numerical Stochastic 
Perturbation Theory (NSPT) JTIJ. We point out that while Xq has been known for a long 
time, X\ has only recently been computed ||11|| , while the only results available for X2 are at 
the moment the ones presented at the Lattice 99 conference [12, 13|, which we will compare 
our result to. We stress that the present computation is still a partial result, since we stay 
within the quenched approximation. We emphasize that an unquenched, non-perturbative 
computation of the binding energy is already available, so that an extension of our result to 
include the fermionic loops contribution is compelling. We are actually working on it. With 
this respect the present computation is a first step which nevertheless we think is useful for 
at least a couple of reasons. First of all, we will show that NSPT can get the result to an 
accuracy which is good enough for the purpose at hand. On top of that we will discuss 
how our result has an impact on the final accuracy on m^m^), at least in the quenched 
approximation. In the next section we briefly summarize the strategy of our computation 
and discuss the numerical analysis which needs to be performed. In section 3. we report our 
results and their impact. Finally in section 4. we state our conclusions. 



2 The computational strategy 

An obvious choice in order to compute the residual mass would be the computation of the 
heavy quark propagator. This would ask for gauge fixing. While our computational tool 
(NSPT) can manage such a task, it turns out that gauge invariant computations are in NSPT 
both easier and less time demanding (for further technical details on this and other points 
concerning the method itself, see |H|). The fact that a gauge invariant computation can 
be set up in order to obtain the residual mass has already been pointed out and exploited 
in ||TTJ and also in [|1^, [Tj|. In order to understand, consider a gauge invariant loop such 
as a Wilson (W) or a Polyakov (P) loop. There are three sources of divergences for such a 
quantity ||I4[. A first one is the logarithmic divergence that can be as usual absorbed in the 
definition of the renormalised coupling. An additional logarithmic divergence appears if the 
contour has an angle (i.e. such a problem is there in the case of a Wilson loop and not for a 
Wilson line or a Polyakov loop); this is usually referred to as the corner divergence. A third, 
linear divergence is the one we are interested in, resulting in (consider for example a Wilson 
loop and note that the divergence exponentiates and shows up as a factor) 

(W) = exp(-cL/a)W lo9 (11) 

In the above formula L is the length of the loop and a the cutoff scale at small distance 
(we have of course already in mind the lattice spacing). The notation Wi og reminds that 
only logarithmic divergences are left. The Wilson line (also called P-line) one would need to 
compute in order to obtain the lattice heavy quark propagator would of course show the same 
divergences. It is interesting to note that HQET was not there in the early eighties when 
people first got interested in the so-called "phase factors" in Gauge Theories. Still of course 
the physical picture was already there, the linearly divergent factor being regarded as mass 
renormalisation of a test particle if one treats the loop as an effective amplitude for a test 
particle moving along the trajectory given by the loop itself. Martinelli and Sachrajda |l 1 
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exploited the above considerations in order to compute 5m to second order by extracting the 
linearly divergent contributions from the whole set of integrals contributing to the Feynmann 
graphs for Wilson loops W(R, T) of generic size R x T in lattice perturbation theory The 
latter ones had already been computed to order in JTJ|]. In |0| our group pinned down 
a first number for 5m to third order (basically its order of magnitude) by performing a fit 
of squared Wilson loops of different sizes L to the expression given in Eq. (|ll|)- At the 
same Lattice 99 conference the Fermilab group presented a value for the same quantity by 
considering a Polyakov loop P in which the mass term contribution was isolated by the 
Coulomb (self-)interaction of the heavy quark Let us now introduce the method we 

exploited in this work. We think that a very neat way to extract 5m goes along the lines we 
now proceed to describe. Consider the quantity 

VW^tmVrW; V T (R) = log ( "ffij^ ) • (12) 

W(R,T) are Wilson loops of size Rx T. Note the dependence on T which is removed by 
taking the T —>■ oo limit; we do not further comment on this point at the moment (there will 
be a lot to say about that later). Eq. (fl2| ) goes of course back to the old subject of Creutz's 
ratios defined to extract the static inter-quark potential. It is for example considered in 
T5fl . It actually achieves the separation of the different (divergent) contributions we are 



interested in. First of all, note that the corner divergences are cancelled between numerator 
and denominator, the number of corners being the same. We now compute everything in 
lattice perturbation theory. In the end one is left with 

V(R)=2 5m + V C oui(R). (13) 

In the above expression one reads the sum of 2 5m and of the actual Coulomb static inter- 
quark potential Vc ou i{R)- The logarithmic divergence connected to the renormalisation of 
the charge is contained in the potential. The latter actually defines a renormalised coupling 
in the so-called potential scheme according to 

Vcom(R) = -C F ^-. (14) 

Note that till now we have made no use of the fact that we are going to compute all the 
quantities in the lattice scheme, apart from having explicitly written the linear divergence 
as 5m. By computing everything in lattice perturbation theory and by separating the linear 
from the logarithmic divergence one actually computes the coupling in the potential scheme 
(av) as an expansion in the lattice coupling (ao), i-e. one computes the matching between 
the couplings in the potential and the lattice scheme. In general at any given order such a 
matching is known once one knows the ratios between the A parameters and the differences 
between non-universal coefficients of the /3-functions {bi,i > 2} in the two schemes. If we 
stop at the order we are interested in 

a 2 {sfi) = ai(jii) + Ci(s) cti(/i) 2 + C 2 (s) a^fi) 3 + . . . (15) 

where 

A 2 



Ci(s) = -26 logs + 26 log- 

Ai 



A 2 b^-bP 



C 2 {s) = C^sf- 2& 1 logs + 2& 1 log-^+ 2 - 2 

Ai Oo 
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Now the point is that at this order the matchings of both the potential |T6| and the lattice 
fT7|| schemes to the MS scheme are known. Since the MS /3-function coefficients are known 
to three loops, one obtains from the above matchings both the ratios of A parameters and 
the 62 coefficients for both the potential and the lattice schemes. This in turns means that 
the matching between the couplings in the potential and the lattice scheme is known to the 
order we need. This defines our strategy. From what we have just said it follows that in 
Eq. ( |13"D everything but 5m is known, i.e. 

V{R) = 2 5m - (a + a(R) a 2 , + c 2 (R) a 3 + ...) (16) 

where 

d(R) = 26 log.R + 26 log-^, log^ = 4.4076 

A h {v) - h i0) h iv) - 
c 2 (R) = c 1 (R) 2 + 2b 1 \ogR + 2b 1 \og^ + b2 - 62 , ' 2 - ' 2 = 8.5794 

A b b 

and 

6m=Y,X^o% ¥1 , (17) 

n>0 

the expression for ay having been taken from Eq. (|l^). The "perfect" strategy would then 
be as follows. Compute V(R) in lattice perturbation theory and compare order by order the 
result to Eq. fli"6|) : the coefficients in the expansion of 5m are the only unknown quantities, 
that could be easily obtained. In practice we had to go through a number of approximations: 
we will discuss how they have to be taken under control and their impact on the final result 
estimated as a final error. First of all, we worked not in infinite volume, but on a lattice of 
fixed 32 4 volume, as it is obviously requested by our computational tool (NSPT). We will 
further comment on this point later. Of course, this implies that we do not have access to 
the real V(R), but only to its (finite volume) Vr(R) approximations. Of course for small R 
and T and for finite lattice spacing the expression in Eq. (|T6|) is distorted by lattice artifacts. 
On top of all the above considerations one should also keep in mind that the results of NSPT 
are affected by statistical errors and their impact should be carefully assessed as well. In 
the end, what we did was to fit the Vr(R) to the asymptotic expression for V(R). The 
approximations Vr(R) were in turn obtained from the computation of W(R,T) for a wide 
range of R and T. A first, obvious point to make is that taking lattice artifacts under control 
asks for having both R and T/R fairly large. The choice of a 32 4 lattice was motivated by 
the compromise between these observations and the computing requirements of NSPT in 
terms of both memory and flops. Having said that, it is also important to point out from 
the very beginning that there are a couple of attitudes one can take with respect to the fits 
we had to perform. First, one can take the coefficients in the expansion of 5m as the only 
unknown parameters, taking for all the other parameters entering Eq. ([16|) their asymptotic, 
lattice artifacts free values. On the other side, one can study the stability of results by fitting 
the other constants in Eq. (|16D as well. This can give an idea of the extent to which our 
computations are good approximations of the infinite volume, lattice artifacts free results. 



3 Results 
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3.1 5m to third order 



We computed W(R,T) on a 32 4 lattice for all R,T < 16. Since we needed to compute 
Sm to order ccq, these expansions had to be computed to the same order. We refer people 
interested in these results themselves to a future review publication |T(J in which they will 



be made available. In that context a detailed discussion of the statistical errors associated 
to NSPT results will be presented. For the purposes at hand here we only need to observe 
that a percent precisionP] of order < 10~ 3 is easy to obtain in the range of R and T we 
were interested in. To do better than this would have been relatively easy, which thing we 
decided not to do in this exploratory study in the quenched approximation. Of course the 
statistical error is magnified by the extraction of the potential, since the log series are there. 
In the end, a percent precision only of order 10~ 2 was obtained for the potential at order 
a.Q, again in the range of R and T we took into account. As far as the choice of the latter is 
concerned, we have already said that reducing lattice artifacts and keeping finite size effects 
under control were the major issues. Note that while the overall volume was kept fixed, both 
the R/T and R scales were relevant. It is trivial to observe that there is no finite size scaling 
that could be set up within these constraints. What we did was to constrain the values for 
R and T this way: for every choice of T in the range 12 < T < 15 the fits were done within 
an interval of 3 -j- 4 values of R subjected to the constraints that a) the central value should 
be such that R/T stayed around a value R/T ~ 3 and b) R > 3. It turns out that with 
such a procedure one is in a good position to perform the fit in which every parameter but 
the coefficients in the expansion of Sm is taken as given by the asymptotic expressions. The 
results so obtained were regarded as approximations of the infinite volume, lattice artifacts 
free results. Of course on top of the statistical errors a systematic error has to be taken 
into account connected to the choice of the fitting interval, i.e. of R and T. Note that as 
approximations of the infinite volume results, the Vr(R) should be loosely dependent on T. 
These can be tested also by studying the stability of results with respect to averaging as 
in (Vt(R) + Vt+x(R))/2. While this gives no new information at first and second order, it 
turns out that it lets you gain something at afj order, where the impact of statistical errors 
(which are bigger at this order) is a bit relieved by this averaging procedure. In order to 
gain some insight into the effectiveness of the procedure we have just described (that is, the 
one in which only Sm is to be fitted), we first of all tried to reproduce the coefficients that 
were already known in the expansion of Sm. Expressing them in an expansion in the bare 
lattice coupling a>o, these coefficients are X = 2.1173, X\ = 11.1520. We got 

JQ = 2.118(2) X[= 11.14(3). 

T The agreement with the analytical, infinite volume results is surprisingly good. The errors 
are in practice dominated by the systematic uncertainty connected to the choice of the fitting 
interval. Such an uncertainty was estimated to be given by the standard deviation pertaining 
to the sample of results of fitting procedures within the constraints we described a few lines 
above. For an extra peace of mind these systematic errors and the statistical ones were added 
not in quadrature. Before making further comments, one can proceed straight to repeat the 
same procedure at order ozq, which gives the result we are interested in as 



X 2 = 86.2(6). 
2 We quote the worst precision, that is the one on order eng. 
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Figure 1: A typical result of the fitting procedure described in the text. The three plots refer 
to first, second and third order of the potential. The solid lines are the appropriate orders 
as given by Eq. (|T6| ) once 8m has been fitted, while circles are the approximations Vt(B). 
In this case T = 14. 



This time an appreciable contribution is also coming from statistical errors in the computa- 
tion of the loops and of the potential, even if the systematic uncertainty connected to the 
choice of the fitting interval nevertheless stays as the major one. Again, these different errors 
were added not in quadrature. In the figure we plot an example of the results one gets by 
this procedure: by fitting values for the coefficients in the expansion of 5m one obtains that 
order by order the approximation of V(R) falls on the curve given by Eq. flltf ). One should 
keep in mind that the effectiveness of the procedure is actually magnified by the figure. This 
in particular true for the highest order, for which the fitted result is a huge shift with respect 
to the final order of magnitude of plotted data. The sensitivity with respect to such a huge 
shift can be easily underestimated if one only looks at the figure. Note by the way that the 
huge shift gives a good idea of the impact of the renormalon on the expansion. Now the 
question comes whether the error bars obtained in the first procedure we described are to 
be fully trusted. As we have already pointed out the effectiveness in reproducing the known 
coefficients within fairly tiny errors is impressive. As it is obvious, this is nevertheless only 
known a posteriori and in this sense it is not a proof, but only a hint at the reliability of the 
method. As stressed, our first procedure relied on the assumption that the Vr(R) were good 
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approximations of the infinite volume, lattice artifacts free V(R), the only sources of errors 
being the statistical ones and the systematic effects of the choice of the fitting intervals. As 
we have already suggested, a second attitude can be taken with respect to fitting the Vr(R) 
to the expression given by Eq. flip]). This does not assume the other coefficients entering the 
formulas as given by their asymptotic, lattice artifacts free values. This means at leading 
order to fit also the coefficient Cf, at second order also the ratio of the A parameters and at 
third order the difference in the /3-function coefficients h 2 as well. This is an indirect test of 
how asymptotic our approximations are. One should anyway always keep in mind that this 
is not a good starting point for a proper finite size analysis, since in all our procedure there 
is a finite fixed volume and a couple of scales around. We will anyway assume the typical 
variation of coefficients with respect to this procedure as a sensible bound on finite size 
and lattice artifacts effects. The constraints we imposed before on the choice of the fitting 
intervals had to be somehow released in order to be able to fit more than one parameter. 
Anyway, no more than two parameters were fitted simultaneously (that is, at every order 
the relevant coefficient in the expansion of 5m plus another constant). This was obtained 
by letting every order inherit fixed values for other parameters. For example: the second 
order fit was in charge of fitting X\ and the ratio of A parameters, inheriting the value for 
Cp either from first order fit or from the infinite volume result (of course also variations 
with respect to these variants were taken into account). Depending again on the choice of 
fitting intervals, we found a rough bound of 10 4- 15% on the variations of all the parame- 

ters entering the Coulomb potential, including 2 bo 2 . Ironically, the latter used to be for 
us a first source of confidence in our results. Given our ignorance of reference [IB|, at the 



beginning we fitted this coefficient together with the third order 6m, assuming for the other 
parameters their asymptotic values. The reasonable agreement of the result of the fit with 
|16|| turned out to be a reasonable hint at the reliability of our result for X 2 . As far as the 



values of the coefficients of the expansion of 5m we found instead a rough bound of 1%. Note 
how this bound appears surprisingly huge with respect to the deviation of our results from 
the exact ones for first and second order. As we have already observed, this is something 
which can be anyway verified a posteriori. A possible step forward in the estimate of finite 
size effects could be to repeat the whole procedure on a different volume, which thing we 
propose to do in the future unquenched study. At the moment, taking into account also the 
last considerations, for an extra piece of mind we assume as our final result 

T 2 = 86.2(0.6) (1.0). 



Our result is quite apart from the one [|13[] presented at the Lattice 99 conference by the 
Fermilab group, suggesting a fairly big impact of systematics. We point out that the com- 
putational tools matter of fact quite different. What we plan to do in the future (for 
the unquenched computation) is to follow also their approach based on Polyakov instead of 
Wilson loops, which thing allows a fine finite size analysis. 



3.2 The impact of the result 

It is now the proper time to address the question of the impact of our result. As a matter 
of fact the result presented in this paper has already been known for some time to people 
involved in the field. Lubicz's review on quark masses at the Lattice 2000 conference |18|| has 
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already referred to the impact of the knowledge of X 2 in the determination of the 6-quark 
mass, at least in the quenched approximation. Basically, things go as follows [I9| (this is 
what one would refer to as a Next to Next to Next to Leading Order - NNNLO - quenched 
calculation of the 6-quark mass). The lines along which the computation can go have already 
been widely presented in the Introduction, the final formula being Eq. (|H]). In order to 
understand what one achieves with the result for X 2 a good point is the comparison with 
the lower order (NNLO) determination^. Of course all the following considerations about the 
impact of NNNLO hold for the quenched case, given the limitation of our result. The residual 
error due to unknown perturbative corrections in Eq. fllPf ) is halved with respect to NNLO 
and the new order improves the convergence of the series. Given that the determination 
of the quark mass can be got from results obtained at various values of (3, a remarkable 
point is that the dependence of the final result for the quark mass on the lattice spacing 
is smaller than in NNLO and very weak. In the end, it turns out that the overall error on 
the mass is of order 1 2%, the contributions being equally split between perturbative and 
non-perturbative ones (this final error was at NNLO dominated by the ignorance of X 2 ). 



4 Conclusions 



We computed the perturbative expansion of the residual mass term in lattice Heavy Quark 
Effective Theory to order oiq. This is an important building block in a renormalon safe 
determination of the 6-quark mass from lattice simulations. Errors appear to be under a 
fairly good control and the impact of the result on the (unquenched) determination of the 
quark mass is quite important. Nevertheless, the result we presented is still a partial one: 
going for unquenching is compelling. The latter is our proposed next step forward. We will 
do that trying to still improve on systematics. 
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